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ABSTRACT 


The  Fast  Field  Program  (FFP)  is  a  technique  for  applying  Fast  Fourier  Transform  methods  to 
Field  theory.  Interim  results  are  presented  for  the  FFP  as  it  v-ould  be  applied  to  a  multilayered 
ocean.  Earlier  work  on  this  subject  has  demonstrated  the  accuracy  of  the  approach,  but  a  reduc¬ 
tion  in  computing  time  was  needed  to  attain  "real  time,”  in  situ  propagation  estimations.  A 
step  in  this  direction  has  been  made  possible  by  the  assumption  that  the  vater  column  is  sub¬ 
divided  into  strata  within  which  the  velocity  of  sound  is  allowed  to  vary  exponentially  with 
depth.  The  kernel  of  the  integral  expression  for  the  pressure  field  is  then  given  in  terms  of  prod¬ 
ucts  of  cylindrical  functions,  which  can  be  calculated  rapidly  through  the  utilization  of  recur¬ 
rence  relations.  Comparisons  of  the  results  of  the  FFP  and  normal  mode  theory  are  provided  for 
both  a  shallow  and  a  deep  water  example  at  low  acoustic  frequencies.  In  the  latter  case,  con¬ 
sideration  is  given  to  the  influence  of  the  subbottom  structure  upon  propagation  predictions.  The 
technique  described  is  not  the  ultimate  solution  for  real  time  calculations;  however,  it  does  rep¬ 
resent  the  fastest  possible  version  of  generating  FFP  predictions  known  to  the  author  at  this  time. 
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The  Fast  Field  Program  (FFP)  is  a  technique  for  applying  Fast  Fourier  Transform 
methods  to  Field  theory.  Interim  results  are  presented  for  the  FFP  as  it  would  be  applied 
to  a  multilayered  ocean.  Earlier  work  on  this  subject  has  demonstrated  the  accuracy  of  the 
approach,  but  a  reduction  in  computing  time  was  needed  to  attain  "real  time, "  in  situ 
propagation  estimations.  A  step  in  this  direction  has  been  made  possible  by  the  assumption 
that  the  water  column  is  subdivided  into  strata  within  which  the  velocity  of  sound  is  allowed 
to  vary  exponentially  with  depth.  The  kernel  of  the  integral  expression  for  the  pressure  field 
is  then  given  in  terms  of  products  of  cylindrical  functions,  which  can  be  calculated  rapidly 
through  the  utilization  of  recurrence  relations.  Comparisons  of  the  results  of  the  FFP  and 
normal  mode  theory  are  provided  for  both  a  shallow  and  a  deep  water  example  at  low  acousti  2 
frequencies.  In  the  latter  case,  consideration  is  given  to  the  influence  of  the  subbottom  struc¬ 
ture  upon  propagation  predictions.  The  technique  described  is  not  the  ultimate  solution  for 
real  time  calculations;  however,  it  does  represent  the  fastest  possible  version  of  generating 
FFP  predictions  known  to  the  author  at  this  time. 
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FAST  FIELD  PROGRAM  FOR  MULTILAYERED  MEDIA 


INTRODUCTION 

The  capability  to  provide  "real  time, "  in  situ  propagation  predictions  is 
highly  desirable  in  many  diverse  areas  of  underwater  acoustics.  It  is  precisely 
this  wide  diversity,  however,  that  makes  it  difficult  to  define  what  real  time  is. 
The  philosophy  employed  here  is  to  develop  the  fastest  possible  method  of  pro¬ 
viding  estimations  and  let  the  eventual  user  make  the  decision  regarding  time 
for  his  particular  application, 

1  2 

The  previous  work  *  on  the  Fast  Field  Program  (FFP),  a  technique  of  ap¬ 
plying  Fast  Fourier  Transform  (FFT)  methods  to  Field  theory,  has  clearly  dem¬ 
onstrated  the  soundness  of  the  concept  and  the  accuracy  of  the  predictions.  In 
thi3  earlier  work,  two  schemes  were  investigated  for  generating  the  input  data 
to  be  transformed.  The  first  involved  a  numerical  integration  of  the  depth  de¬ 
pendent  differential  equation  for  an  arbitrary  sound  velocity  profile.  In  the 
second  approach,  the  variation  of  the  sound  velocity  with  depth  was  restricted 
in  such  a  way  that  known  solutions  were  available.  Intrinsic  difficulties  of  the 
same  type  were  found  for  both  methods,  and  it  was  felt  that  the  required  com¬ 
puting  time  could  be  reduced.  The  technique  described  in  this  report  offers 
such  a  reduction. 


THE  FAST  FIELD  PROGRAM  (FFP)  AND  THE 
FAST  FOURIER  TRANSFORM  (FFT) 

o 

It  is  well  known  that  the  field  due  to  a  monochromatic  point  source  can  be 
represented  as 

r°° 

*(z ,  r)  =  G(z ,  zs;  *)  H ^  «r)  *  d$  .  (1) 

—OD 

The  scalar,  4*  ,  is  related  to  the  pressure  field  and  the  particle  velocity  through 
the  expressions 


where  p  is  an  abbreviation  for  the  density.  The  function  G  represents  the 
depth  dependent  Green's  function,  which  must  simultaneously  satisfy 


PS’S®*! 


and  the  input  to  the  FFT  is  found  to  be 


t?  _  rt*  ,  .  t  \  imroAJ 

E  -  u(z,z  ;  £  )  £  e 
m  s  m  m 
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Equation  (4)  is  now  in  a  form  that  is  directly  amenable  to  the  FFT.  The 
result  of  the  application  of  the  FFT  is  the  value  of  the  field  at  each  of  the  dis¬ 
crete  N  ranges.  Once  the  input,  Eq.  (5),  is  obtained,  the  calculation  of  Eq. 
(4)  proceeds  very  rapidly.  For  the  UNIVAC  1108  and  N  =  8192,  this  calculation 
consumes  2.689  sec.  The  main  concern  of  this  report  is  with  the  rapid  calcu¬ 
lation  of  accurate  values  of  the  input.  The  basis  of  a  scheme  that  allows  this  to 
be  done  is  discussed  in  the  next  section. 


STRUCTURING  OF  THE  KERNEL  FOR  THE  RAPID 
CALCULATION  OF  THE  INPUT 

The  wave  number,  k  =  >/(«/cj)&  + 2ia,  in  Eq.  (2)  is  assumed  to  be  a  com¬ 
plex  with  attenuation  coefficient  a  >  0  since  the  Green’s  functions  could  be 
singular  along  the  axis  of  integration.  The  homogeneous  form  of  Eq.  (2)  is  then 
written  as 

<6» 


The  subscript  j  is  included  in  anticipation  of  the  subdivision  of  the  velocity 
profile  into  strata.  In  order  to  provide  insight  into  the  motivation  for  this  sub¬ 
division,  a  slight  digression  is  necessary. 

The  general  solution  of  Eq.  (6)  for  an  arbitrary  sound  velocity  profile 
can  be  obtained  by  numerical  integration.  Two  linearly  independent  solutions 
can  be  obtained  for  a  given  value  of  w  these  solutions  satisfy  the  relevant 
boundary  conditions  and  thus  allow  the  construction  of  the  Green’s  function. 
These  solutions  change  their  character  at  the  turning  point.  The  transition 
from  an  oscillatory  to  an  exponential  behavior,  or  vice  versa,  can  be  extremely 
sensitive  to  round-off  errors  for  certain  profiles.  At  the  same  time,  the  solu¬ 
tions  exhibiting  a  turning  point  as  a  function  of  $m  with  z  fixed  are  the  most 
important  contributors  to  the  integration  in  Eq,  (1).  It  is  certainly  true  that  the 
accuracy  of  the  process  may  be  improved  by  resorting  to  double  precision  or 
special  integration  techniques.  The  price  paid  for  these  manipulations,  how- 
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ever,  is  additional  computer  time.  Even  if  one  were  to  ignore  this  aspect  of  the 
problem,  the  prospect  still  remains  that  the  integration  must  bo  repeated  for 
each  of  the  N  discrete  values  of  The  result  is  that  the  amount  of  compu¬ 
ter  time  needed  to  generate  the  input  is  unsatisfactory  when  compared  with  the 
time  required  to  do  the  FFT. 

If  the  variation  of  the  profile  is  restricted,  known  solutions  of  Eq.  (6)  are 
available.  This  approach  was  utilized^  by  fixing  the  value  of  the  profile  at  three 
arbitrary  points,  and  the  solutions  were  given  in  terms  cf  hypergeometric  func¬ 
tions.  This  technique  has  essentially  the  same  limitations  as  numerical  inte¬ 
gration.  As  the  frequency  or  water  depth  is  increased,  there  is  a  correspond¬ 
ing  increase  in  the  required  computing  time  and  in  the  difficulty  with  round-off 
errors.  As  was  mentioned  earlier,  when  speed  is  not  the  prime  concern  the 
problem  of  round-off  errors  can  usually  be  circumvented,  and  thus  the  above- 
mentioned  approaches  are  possible  candidates  for  the  generation  of  the  input. 

To  return  to  the  technique  of  this  report,  let  the  profile  be  restricted  to 
vary  exponentially  with  depth  within  each  layer  according  to 

°j(zi,=ci-i<zi)e±<ZrZi)/Hi  zj  -  zi  -  ' 

where  j  =  1,2,  •  •  • ,  and  cQ(z1)  is  the  surface  velocity.  The  z  axis  is  taken  to 
be  positive  in  the  downward  direction.  The  solutions  of  Eq.  (6)  are  then  found 
to  be 

w-c,  blM-  (7) 

j  J 

that  is,  in  terms  of  cylindrical  functions  of  complex  order 

v.  S(£H.  -i«H.)  (8) 

J  J  J 


and  real  argument 


wH. 


e±(zi-zj)/Hj 


(9) 


The  possible  values  of  z.  are  either  one  of  the  interface  depths  or  the  source 
or  receiver  depths. 
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The  variation  of  G  as  a  function  of  £  may  now  be  obtained  economically 
and  rapidly  through  the  recurrence  relations  for  products  of  cylindrical  func¬ 
tion.  This  aspect  of  the  technique  will  be  given  more  detailed  attention  in  the 
next  section.  The  fact  that  the  input  can  now  be  obtained  by  utilization  of  a  re¬ 
currence  relation  is  most  desirable*  but  it  does  not  complete  the  requirements. 
The  recurrence  relation  must  also  be  linear  in  the  variable  £  since  the  sam¬ 
ples  in  the  FFT  need  be  equispaced.  One  other  profile  that  meets  these  speci¬ 
fications  is  known  to  the  author.  Consider  the  equation 


The  basis  of  solutions  is  found  to  be 


/31(z)=P^  jtanh(  JT^)] 

^2(z)  =  Qj‘[tanh(— jj-i)]  ,  (10) 

and  the  desired  recurrence  relation,4  satisfied  by  both  of  the  solutions  in  Eq. 
(10),  is 

u42  x  u+1  u 

p;  (X)  +2(/i+  1)-  - - (x)  +  (1/  —  M)  O'  +  M  +  1)  P„  (x)  =  0  . 

\1  -  X 


The  cosh  profile  of  the  associated  Legendre  functions,  Eq.  (10),  is  not  signifi¬ 
cantly  different  from  the  exponential  profile;  thus,  it  was  decided  to  use  the 
cylindrical  functions. 


THE  SINGLE  EXPONENTIAL  LAYER 

The  case  that  requires  the  least  amount  of  computer  time  is  that  in  which 
the  velocity  profile  can  be  represented  by  a  single  exponential  layer,  as  de¬ 
picted*  in  Fig.  1.  The  bottom  was  assumed  to  be  a  homogeneous  liquid,  semi¬ 
infinite  in  extent,  and  with  a  velocity  faster  than  that  of  the  overlying  water. 

The  densities  of  the  water  and  bottom  were  taken  to  be  identical  for  convenience. 


c(»  =  c0*’z/H  (ft/sec) 


4982  5000 


4>(z)  =  cv(u),  yf  -  k0  ,  v  s  £h  ~  iaH 


Fig.  1.  Environmental  Description  for  a 
Single  Exponential  Layer 
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It  is  faster,  from  a  computational  point  of  view,  to  allow  the  velocity  in  the 
bottom  to  also  vary  exponentially  with  depth.  The  information  could  then  be  ob¬ 
tained  from  a  recurrence  relation  rather  than  the  evaluation  of  a  square  root  of 
a  complex  number  and  a  trigonometric  function  for  each  of  the  N  values  of  £m. 
The  normal  mode  program  that  was  employed  for  comparison  purposes,  how¬ 
ever,  could  not  accommodate  such  a  structure. 


THE  GREEN’S  FUNCTION 

The  constituents  of  the  Green's  function  are  found  to  be 


(fb(YV 


G<VV{)" 


[Q,(7®,^)  -  Y22  P„  (7®, 7*)] 
[Qt,(7j,721)-Y^Pl/(7j,7^)] 


[Q„(7*,72)  -y22  p„  (7^,  y2)] 

[Ql)(7^,72)  -  Y2  P,(7j,72)] 


z„  <  z  <  z 
1 _  r~  s 


z  <  z  <  z  , 
s  r  -  2 


where  the  products  of  cylindrical  functions  are  defined  as 

P„  (a,  b)  =  Jv  (a)  Y„  (b)  -  J„  (b)  Y„  (a) 

Q„(a,b)  =  J„(a)  Y'(b)  -  J' (b)  Y„  (a) 

Rt.  (a,b)  =  J'v  (a)  Yv  <b)  -  J,  (b)  Y'  (a) 

Sj,  (a,  b)  =  J’  (a)  Yl  (b)  -  (b)  Y'  (a)  .  (11) 


The  primes  will  always  denote  differentiation  with  respect  to  the  argument. 

The  last  two  expressions  of  Eq.  (11)  are  not  needed  explicitly  for  this  example 
but  will  be  required  in  the  next  section.  The  admittance  of  the  bottom,  y|  ,  is 
found  to  be 


Ul  <  w/°B 
Kl  >W/CB  • 
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The  representation  for  the  Green's  function  in  terms  of  the  products  of  cy¬ 
lindrical  functions,  Eq.  (11),  is  not  only  symbolically  efficient  but  it  also  rep¬ 
resents  a  savings  in  computing  time.  These  products  can  be  evaluated  rapidly 
through  the  recurrence  relations,  ^ 


Q 

It 

P 

S 


v+1 

y+1 

v+2 

•>+1 


-  R„  +  (u/a)  P„ 


(j£ l)p 
b  v+1 


-  Q„  +  (u/b)  P„ 


(*+1) 

a 


P 


*+l 


P„  -  2  (.+1)  [Q„+1/a  +  R +1/b] 

P.+2  +  (‘'+1»[Q.+l/a+R.+lA]- 


ta!P 

ab  v+1  ' 


(12) 


Two  areas  of  concern  associated  with  Eq.  (12)  are  the  availability  of  start¬ 
ing  values  and  the  stability  of  the  recurrence  process.  For  the  examples 
■examined,  the  WKB  approximations  of  the  functions  defined  in  Eq.  (11)  pro¬ 
vide  accurate  starting  values.  The  question  of  stability  is  most  important  when 
the  value  of  the  real  part  of  v  lies  between  the  arguments  a  and  b.  The  dif¬ 
ference  between  these  arguments  increases  with  either  incre  asing  frequency 
or  water  depth.  The  behavior  of  the  functions  in  Eq.  (11)  within  this  transition 
region  takes  the  form  of  an  oscillatory  function  modulated  by  a  growing  expon¬ 
ential  envelope.  As  this  separation  distance  increases,  the  envelope  grows 
proportionally  and  thereby  enhances  the  possibility  that  round-off  error  will 
affect  the  calculations  of  Eq.  (12).  A  more  detailed  account  of  these  aspects 
and  a  method  of  overcoming  the  latter  is  provided  in  Reference  8. 


SAMPLING  THE  INPUT 

The  sampling  distance  in  the  wave  number  domain  is  determined  by  the 
scaling  factor  H  and  the  fact  that  the  change  in  v  must  be  unity.  Then  from 
Eq.  (8)  one  sees  that 


(13 

The  magnitude  of  H  for  most  applications  will  be  on  the  order  of  10^.  The 
sampling  distance,  Eq.  (13),  may  be  too  coarse  at  some  frequencies;  however, 


8 


it  may  be  altered  as  will  be  explained  in  the  section  on  the  two-layer  model . 

The  parameter  a  is  assigned  the  value  o£  the  attenuation  coefficient  at  the  fre¬ 
quency  of  interest.  The  sampling  region  should  encompass  that  portion  of  the 
real  £  axis  near  the  region  in  which  the  Green’s  function  is  singular.  If  kmax 
is  the  wave  number  corresponding  to  the  minimum  velocity  of  the  water-bottom 
combinat^Ti,  the  beginning  of  the  sampling  region  is  found  from 

L  =  k  -  NA£  .  (14) 

0  max 


The  integrand  decays  to  zero  for  values  of  £  greater  than  kmax .  The  selec¬ 
tion  of  N  then  determines  the  length  of  the  axis  sampled.  The  question  of  how 
large  N  must  be  is  the  analog  of  how  many  modes  one  includes  in  normal  mode 
theory.  It  is  worth  emphasizing  that  the  FFP  is  not  an  eigenvalue  problem.  It 
is  not  necessary  to  know  the  exact  location  of  the  singularities  in  the  complex  £ 
plane,  but  rather  only  the  bounds  of  the  region  where  they  are  distributed. 
Certain  crude  but  useful  techniques  are  available  to  answer  questions  regarding 
&£  and  N.  The  value  of  &£  can  be  doubled  and  N  reduced  to  one-half  of  its 
previous  value.  In  this  way  the  length  of  the  region  sampled,  NA£ ,  is  un¬ 
changed  and  the  effect  of  a  smaller  sampling  distance  can  be  determined.  A 
doubling  of  N  with  &£  unchanged  will  show  the  significance  of  effectively 
adding  more  modes. 


COMPARISON  OF  THE  FFP  AND  NORMAL  MODE  THEORY 

One  method  to  demonstrate  the  relationship  between  the  FFP  and  normal 
mode  theory  is  to  examine  the  magnitude  of  the  input  as  a  function  of  the  wave 
number.  This  information  is  presented  in  Fig.  2  for  the  problem  under  con¬ 
sideration.  The  series  of  peaks  is  the  result  of  nearby  singularities  of  the 
Green’s  function,  which  would  give  rise  to  the  normal  modes.  The  numbers 
associated  with  these  peaks  are  the  mode  numbers  derived  by  Bartberger  and 
Ackler  by  using  their  three-layer,  normal  mode  program.9  In  all,  they  calcu¬ 
lated  37  modes  with  their  last  eigenvalue  having  a  real  part  of  1.217.  Excellent 
agreement  was  found  between  the  location  of  the  peaks  and  the  real  part  of  their 
eigenvalues.  As  can  be  seen,  the  region  sampled  by  the  FFP  includes  more  of 
the  high-order  modes  but  omits  the  first  four  low-order  modes. 

Propagation  loss  versus  range  is  displayed  in  Fig.  3.  The  upper  graph 
represents  the  normal  mode  predictions  and  the  lower  graph  represents  that  of 
the  FFP.  The  range  resolution  for  the  FFP  was  four  times  as  small  as  that 
used  in  the  normal  mode  program,  and  this  fact  accounts  for  the  greater  detail 
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MAGNITUDE  OF  INTEGRAND  (X10 


016 


WAVE  NUMBER  (X 10  ) 


12.55  12.58 


Fig.  2.  Magnitude  of  Input  for  the 
Single  Exponential  Layer 
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Fig.  3.  Comparison  of  FFP  and  Normal  Mode  Predictions:  Single  Layer 


in  that  pattern.  At  those  points  where  the  ranges  of  both  solutions  coincide,  the 
difference  between  the  solutions  is  no  greater  than  1  dB.  This  is  the  reason  for 
plotting  the  two  results  in  a  displaced  manner.  The  FFP  provides  results  out  to 
about  200  kyd  in  increments  of  26  yd  with  the  input  values  given  in  Fig.  2.  The 
total  FFP  results  are  shown  in  Fig.  4,  where  the  effect  of  aliasing  in  the  FFT 
calculations  is  evident  at  the  longest  ranges. 

The  total  execution  time  for  this  example  was  31  sec  on  the  UNIVAC  1108. 
The  ratio  of  time  needed  to  generate  the  input  to  the  time  required  to  actually  do 
the  FFT  is  thus  about  nine  to  one .  Of  all  the  techniques  considered  for  the  gen¬ 
eration  of  the  input,  ;he  single  exponential  layer  requires  the  least  amount  of 
computing  time  and  is  considered  to  be  optimum  at  this  stage  in  the  develop¬ 
ment. 

Unfortunately,  a  single  exponential  layer  is  not  a  representative  environ¬ 
mental  model.  The  possibility  of  using  multilayers  does  exist,  however,  and 
the  ramifications  of  this  approach  are  discussed  in  the  next  section. 
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THE  TWO-LAYER  MODEL 


A  typical  warm -weather  profile  for  the  Tyrrhenian  Sea  is  given  in  Fig.  5. 
The  velocity  of  sound  versus  depth  can  be  represented  by  the  two  exponential 
layers, 


Cjtej)  =  c0(Zj)  e'*Zl'zl)/,Hl  Zj<  Zj  <  z2 

C2<V  -  °i<z2>  e(Zi'Z2,/H2  \  <  *3  • 


The  velocity  at  the  surface,  Cq(z^),  was  taken  to  be  5040  ft/sec.  The  layer 
depths  Z2  and  Z3  are  200  and  9000  ft,  respectively,  and  the  scaling  factors  for 
each  layer  are 


Hx  =  12  x  io3 

H2  =  36  x  io4  . 


Consideration  was  given  to  the  subbottom  structure,  since  the  frequencies 
of  interest  were  in  the  seismic  range.  The  two-layered,  liquid  bottom  finally 
adopted  (Fig.  5)  was  arrived  at  by  trying  to  match  experimental  bottom  loss 
measurements  for  this  area.  The  results  of  those  measurements  are  de¬ 
picted  in  Fig.  6.  The  bottom  loss  in  the  75-  to  150-Hz  band  provided  by  the 
two-layered  bottom  is  given  in  Fig.  7.  It  can  be  seen  that  the  theoretical  curve 
contains  all  the  important  features  present  in  the  measurements. 
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VELOCITY  (ft/sec) 

5010  5050  5090 


Environmental  Description  for  the  Tyrrhenian  Sea 
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INTEGRATION 


THE  GREEN'S  FUNCTION 


If  the  same  abbreviations  introduced  in  the  previous  section  are  utilized, 
the  Green's  function,  G(zs,zr;£),  for  the  case  of  the  source  in  the  first  layer, 
is  found  to  be 


.  1  r 

.W 


z,  <  z  <  z 
1~  r  -  s 


-7rH  1  S 

— ip„  (VV 
2  1 


-(t)Y2 


z  <  z  <  z 
s  ~  r  - 


2 


16 


(15) 


The  admittance  at  each  of  the  layer  interfaces  is 


where 


p  -H, 

r ii 


y2  -(xb - — ;  Tir— — — r  2=22  <10> 


2 


K1  £K2  "  kl tanh  {ikz  d)] 
3  ~  [Ki  ”  k2  tanh  (ik^  d)  ] 


Y?  = 


z  =  z 


3  ’ 


(17) 


$2  +  2ia_, 

K.  -*■*,  '  B  j = i»2  , 

J  (wPj) 


+  2ia, 


B 


d  is  the  thickness  of  the  first  bottom  layer,  and  subscripts  w  and  B  refer  to 
water  and  bottom,  respectively. 
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MULTILAYERS  AND  SAMPLING 


The  addition  of  multilayers  and,  thus,  different  scaling  factors  provides  a 
choice  for  the  fundamental  value  of  ,  the  sampling  distance.  The  facts  that 
the  FFT  requires  equispaced  inputs  and  that  the  value  of  is  different  for 
each  layer  would  seem  to  be  in  conflict.  This  situation  can  be  resolved  by  set¬ 
ting  equal  to  the  reciprocal  of  the  largest  scaling  factor  and  then  taking  the 
other  scaling  factors  to  be  integer  multiples  of  this  value.  In  this  instance,  one 
has 


As  before,  the  beginning  cf  the  sampling  region  is  found  from  Eq.  (14).  The 
discrete  values  of  the  order  for  each  layer  are  then  found  to  be 


<'2)m"«0H2-ka2)  +  m 


(18) 


H1 

^m  <*0H1  ■  +  m  H  * 


(19) 


since 


t  =  $  +JL 

m  0  V 


where  m  =  0,1,2,  •••  ,N-1. 

The  recurrence  relations  (Eq.  (12) )  yield  the  value  of  the  required  func¬ 
tions  at  successive  values  of  the  order  which  differ  by  unity.  Those  values 
would  coincide  exactly  with  the  discrete  values  of  Eq.  (18)  but  not  those  of  Eq. 
(19).  The  remedy  for  the  latter  case  is  to  have  H2/H1  starting  values  for  the 
first  layer.  The  relationship  between  the  variable  of  integration,  the  order  for 
each  layer,  and  the  path  through  these  values  when  utilizing  Eq.  (12)  is  illus¬ 
trated  graphically  in  Fig.  8.  The  implication  of  this  process  on  computer  stor¬ 
age  is  provided  in  Table  1. 
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Table  1 


VARIABLES  AND  REQUIRED  LOCATIONS 


Variable 

Required  Locations  (Complex) 

P 

2 

1 

Q 

^  V 

1 

1 

R 

1 

1 

S 

1 

1 

P 

60 

2 

Q 

30 

2 

R 

30 

2 

S 

30 

2 

G(z  ,z  ;£  ) 

N 

s  r  m 

It  is  possible  to  set  A£  equal  to  the  reciprocal  of  the  smallest  scaling  fac¬ 
tor,  calculate  Eq.  (12)  for  v£  as  before,  but  only  retain  every  thirtieth  value. 
This  procedure  results  in  a  larger  fundamental  sampling  distance.  However, 
this  result  is  also  available  with  the  first  scheme.  Suppose  the  smallest  funda¬ 
mental  sampling  distance  is  decided  upon.  If  it  is  desirable  to  double  this  value, 
one  would  merely  have  to  do  twice  as  many  calculations  and  retain  every  other 
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one.  On  the  other  hand,  the  sampling  distance  can  be  halved  by  repeating  the 
procedure  outlined  in  Table  1  twice  with  N  =  N/2.  The  second  time  the  start¬ 
ing  point  would  be  £0  =  £0  +  A£/2 . 

COMPARISON  OF  THE  FFP  AND  NORMAL  MODE  THEORY 

The  normal  mode  calculations  were  performed  by  Bartberger  and  Ackler  at 
NADC  using  their  three-layer  program.  A  modification  of  the  program  was  re¬ 
quired  to  accommodate  a  multilayered  bottom.  This  was  accomplished  H  by 
retaining  only  the  magnitude  of  the  reflection  coefficient  at  the  water  bottom. 

The  magnitude  of  the  reflection  coefficient  was  that  produced  by  the  two-layered 
bottom  used  in  the  FFP  calculations.  All  normal  mode  predictions  associated 
with  this  example  were  obtained  with  the  above-described  reflection  coefficient. 
FFP  calculations  were  performed  with  the  two-layered  bottom  and  with  a  semi¬ 
infinite  bottom  having  the  velocity  of  the  second  bottom  layer.  By  utilizing  both 
types  of  bottoms  and  increasing  the  frequency  while  other  parameters  were  held 
constant,  it  was  hoped  that  the  significance  of  the  subbottom  structure  could  be 
determined. 

A  qualitative  indication  of  the  influence  of  the  bottom  can  be  obtained  from 
an  examination  of  the  Green's  function.  Consider  Eq.  (15)  for  the  case  when  the 
receiver  is  in  the  second  layer.  The  smallest  argument  of  the  product  functions 
for  the  second  layer  is  72  ,  since  the  velocity  gradient  is  positive  in  this  layer. 
The  function  Y„2(72)  will  eventually  dominate  these  product  functions  when 
ii/2i>7sj .  At  some  point  in  the  sampling  process  one  has  the  approximate  values 


and,  similarly, 


P,  <VT2,-J-  (r2,Y.  <T2> 

2  2  2 


Q,  <vt2»  sj, Yiji> 

La  La  La 


%(T2-Y2»S,,2<Y2,%<Y2) 


<YT2)SJ.  <T2)Y:  (12> 

2  2  2 
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Substitution  of  these  expressions  into  the  Green's  function  produces 


G(z 


P,  <VTi) 

1 


LI 

"  ,  1  2k  /  1\,2  _  1  2* 


v? 


(20) 


where  the  admittance  of  the  second  layer  reduces  to 


2  \-H  / 


2  Ji  <y 
2  2 


2  \  <y 
2  2 


This  expression  is  now  independent  of  y|  ,  and  thus  the  bottom  has  no  in¬ 
fluence  at  these  values  of  ^  *  When  the  receiver  is  in  the  first  layer,  a  simi¬ 
lar  reduction  of  the  Green's  function  is  possible.  As  one  might  expect,  the 
values  of  the  order  for  which  this  occurs  are  related  to  the  rays  that  would  have 
a  bottom  vertex  within  the  water  column. 

A  comparison  of  normal  mode  and  FFP  predictions  is  provided  in  Figs.  9 
through  14  for  various  frequencies.  It  is  worthwhile  to  repeat  that  all  the 
normal  mode  calculations  were  made  with  the  two-layered  bottom  and  that  the 
FFP  results  are  either  for  the  same  two-layered  bottom  or  for  the  semi¬ 
infinite  bottom.  An  indication  of  the  influence  of  the  subbottom  is  most  evident 
in  Fig.  9  at  25  Hz,  the  lowest  frequency  for  which  comparisons  were  attempted. 
In  this  figure,  the  semi -infinite  bottom  was  used  in  the  FFP.  The  agreement 
between  the  two  models  is  good  to  about  40  kyd.  The  peak  vaiues  of  the  normal 
mode  solution  beyond  this  point  show  a  loss  of  about  5  dB  less  than  their  FFP 
counterparts.  In  addition,  it  can  be  seen  that  the  two  beat  interference  pattern* 
are  dissimilar.  Some  of  the  energy  entering  the  two-layered  bottom  would  be 
returned  to  the  receiver  after  reflecting  from  the  second  bottom  interface, 
whereas  it  is  never  returned  with  the  semi-infinite  bottom.  This  hypothesis 
was  tested  by  including  the  two-layered  bottom  in  the  FFP  calculations.  The 
results  of  that  calculation  are  displayed  in  Fig.  10.  The  pattern  for  the  FFP 
at  close  ranges  has  the  same  broad  features  as  those  in  the  previous  figure. 


*  Evidence  of  this  fact  is  provided  in  Fig.  15  and  will  be  discussed  at  a  later  point  in  this 

section. 
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Propagation  Loss  for  a  Semi-infinite  Bottom  and  a 
Frequency  of  25  Hz 
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Fig.  10.  FFP  Propagation  Loss  for  a  Two-layered  Bottom  and  a 

Frequency  of  25  Hz 


However,  a  more  detailed  structure  is  present  and  suggests  the  presence  of 
stronger  bottom  arrivals.  This  perturbation  diminishes  with  increasing  range, 
a  fact  that  could  be  accounted  for  by  the  decibel  loss  per  bounce  analog  in  ray 
theory.  A  further  comparison  of  the  two  figures  reveals  that  the  peak  values  of 
both  results  are  in  better  agreement.  , 

Similar  comparative  runs  were  made  at  35,  50,  and  100  Hz  to  determine 
the  impact  of  the  subbottom  with  increasing  frequency.  As  might  be  expected, 
the  subbottom  becomes  less  significant  as  the  wavelength  decreases.  A  com¬ 
parison  of  the  FFP  and  normal  mode  theory  is  provided  in  Figs.  11  and  12,  at 
35  and  50  Hz,  respectively,  for  the  two-layered  bottom.  The  results  at  100  Hz, 
with  the  different  bottom  types,  are  given  in  Figs.  13  and  14. 

A  comparis  in  of  these  two  figures  reveals  that  at  close  ranges  the  pattern 
for  the  two-layered  bottom  (Fig.  13)  contains  a  finer  structure  than  that  for  the 
single-layered  bottom.  This  effect  is  the  same  as  that  observed  at  25  Hz. 

Since  no  systematic  difference  in  level  exists,  however,  it  appears  that  at  100 
Hz  most  of  the  energy  is  trapped  in  the  water  column. 

Isolated  ranges  exist  at  all  frequencies  where  the  FFP  and  normal  mode 
predictions  are  in  wide  disagreement.  This  was  not  the  case  in  the  single-layer 
comparisons  of  Fig.  3.  It  is  felt  that  these  isolated  discrepancies  are  partially 
due  to  the  omission  of  the  phase  information  of  the  reflection  coefficient  in  the 
normal  mode  calculations.  A  second  difference  in  the  two  results  occurs  at  the 
very  close  ranges  where  the  FFP  results  for  the  two-layered  bottom  show  sub¬ 
stantially  less  loss  than  the  normal  mode  results  for  all  the  frequencies  presented. 
Explanations  in  terms  of  such  things  as  (1)  the  phase  of  the  reflection  coefficient, 
(2)  branch  cut  contributions,  and  (3)  the  fact  that  the  profile  varies  exponentially 
with  depth  for  the  FFP  and  linearly  in  the  normal  mode  program,  exist  but  have 
not  been  substantiated. 

Additional  evidence  of  the  contribution  from  the  layered  bottom  is  provided 
in  Fig.  15,  where  the  magnitude  of  the  input  at  35  Hz  is  plotted  versus  the  var¬ 
iable  of  integration.  The  results  for  the  semi-infinite  and  two-layered  bottom, 
respectively,  are  plotted  on  alternate  levels  for  comparative  purposes.  The 
magnitude  is  larger  for  the  two-layered  bottom  over  the  entire  region  sampled 
except  at  the  very  end.  The  larger  peaks  correspond  to  modes  associated  with 
the  bottom.  The  fact  that  the  two  results  are  almost  identical  for  wave  num¬ 
bers  greater  than  .  043  is  explained  by  means  of  Eq.  (20);  that  is,  in  this  region, 

I  *'2l>'y2  the  Green's  function  is  independent  of  the  bottom  admittance. 
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Fig.  12.  FFP  Propagation  Loss  for  a  Two-layered  Bottom  and  a 

Frequency  of  50  Hz 
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Fig.  13.  FFP  Propagation  Loss  for  a  Two-layered  Bottom  and  a 

Frequency  of  100  Hz 


Propagation  Loss  for  a  Semi-infinite  Bottom  and  a 
Frequency  of  100  Hz 
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Fig.  15.  Comparison  of  Input  at  35  Hz  for  Different  Bottom  Structures 
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Fig.  15  (Cont'd).  Comparison  of  Input  at  35  Hz  for  Different  Bottom  Structures 
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Fig.  15  (Cont’d).  Comparison  of  Input  at  35  Hz  for  Different  Bottom  Structures 


An  8192  complex  point  FFT  was  U3ed  in  all  the  above  cases.  The  sampling 
distance,  ,  was  taken  as  the  reciprocal  of  H2  at  25  and  35  Hz.  At  50  and 
100  Hz  it  was  set  equal  to  twice  this  value  in  order  to  provide  a  larger  sam¬ 
pling  region.  At  the  lower  frequencies,  the  FFP  yields  propagation  loss  esti¬ 
mations  out  to  about  740  kyd  in  intervals  of  92  yd.  It  is  difficult  to  establish  the 
range  at  which  the  effects  of  aliasing  in  the  FFT  are  significant.  A  conserva¬ 
tive  estimate  would  be  the  midpoint,  which  for  these  examples  is  180  nautical 
miles.  For  this  reason,  only  one-half  of  the  results  have  been  presented  in  the 
figures  in  this  section.  The  execution  time  for  the  UNIVAC  1108  was  roughly 
98  sec.  Of  this  total  time,  2.689  sec  is  required  to  do  the  FFT.  Thus,  the 
ratio  of  input  generation  time  to  FFT  time  is  about  30  to  1  for  this  two-layer 
problem. 


DISCUSSION 

Interim  results  on  the  development  of  the  Fast  Field  Program  (FFP)  have 
been  reported.  It  is  felt  that  the  technique  discussed  offers  a  major  reduction 
in  execution  time  over  previously  investigated  approaches.  This  has  beei  made 
possible  by  describing  the  sound  velocity  profile  in  terms  of  an  exponential  var¬ 
iation  with  depth,  which  permits  the  rapid  calculation  of  the  input  to  the  FFT 
through  the  use  of  recurrence  relations  for  products  of  cylindrical  functions. 

The  required  execution  time  in  utilizing  this  technique  for  a  single  expon¬ 
ential  layer  is  the  fastest  possible  of  the  various  methods  of  input  generation 
investigated.  A  more  complicated  environment  may  be  accommodated  via  the 
combination  of  exponential  layers  of  different  gradients.  The  execution  time 
increases  with  each  additional  layer,  but  the  technique  is  still  superior  to  the 
approaches  used  in  the  earlier  development  of  the  FFP.  Along  these  lines  it 
should  be  noted  that  although  the  ratio  of  generation  to  execution  times  was 
thirty  to  one  for  the  two-layer  case,  the  increase  in  time  over  the  single-layer 
example  was  not  entirely  covered  by  the  addition  of  a  second  layer.  The  single¬ 
layer  example  pertained  to  shallow  water,  whereas  the  two-layer  model  was 
for  deep  water.  The  water  depth  has  a  significant  bearing  on  the  total  execution 
time,  since  an  increase  in  water  depth  will  result  in  a  larger  transition  region, 
as  discussed  in  the  section  on  the  Green's  function,  and  thus  will  require  great¬ 
er  computational  accuracy  in  evaluating  the  kernel.  This  point  is  discussed  in 
detail  in  Reference  8. 

Comparisons  of  propagation  loss  predictions  for  the  FFP  and  normal  mode 
theory  show  that  not  only  will  the  FFP  provide  the  results  quickly  but  also  ac¬ 
curately.  A  similar  comparison  of  execution  times  between  these  two  methods 
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of  prediction  could  not  be  made  directly  because  the  results  were  obtained  on 
different  computers. 

In  addition  to  demonstrating  the  FFP  for  a  multilayered  medium,  the  two- 
layer  model  also  serves  to  illustrate  the  need  for  accurate  bottom  loss  informa¬ 
tion  at  low  frequencies.  Inasmuch  as  the  low-frequency  end  of  the  spectrum  is 
attracting  more  attention,  ic  would  be  desirable  to  have  additional  bottom  loss 
information  available. 
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